function [ errors ] = 	EquationError( model_sol )
% Check the accuracy of the solution by evaluating the Euler equation errors.
% model_sol is the post_processing solutions. 

if( ~isfield(model_sol, 'productivity_matrix') )  % check whether the most solution is post processing
    model_sol = post_processing( model_sol );
end

par = model_sol.par;

AH = par.AH;  AL=par.AL;   beta=par.beta;  alpha=par.alpha; 

w = model_sol.w_matrix;
lambda = model_sol.lambda_matrix;

p = model_sol.p_matrix;

xK = model_sol.xK_matrix;
yK = model_sol.yK_matrix;
sigmaK = par.sigmaK;
sigmap = model_sol.sigmap_matrix;

kappaw = model_sol.kappaw_matrix;
kappap = model_sol.kappap_matrix;
kappad = model_sol.kappad_matrix;
kappa_fs = model_sol.kappa_fs_matrix;

delta_x = model_sol.delta_x_matrix; 
indicator = model_sol.indicator_matrix;

muR = model_sol.muR_matrix;
rd = model_sol.rd_matrix;

% check several equations: (1) The first order condition of bankers. (2) The first order condition of households.  (3)
% Equilibrium jumps

% (1) FOC of banks on capital
residuals1 = muR + AH./p - rd - ( sigmaK + sigmap  ).^2.* xK - lambda.* ( (kappap+ indicator.*alpha.*beta)./( 1-xK.*kappap-alpha.*delta_x ) ); 
residual_tes = muR + AH./p - rd - ( sigmaK + sigmap  ).^2.* xK - lambda.* ( (kappap+ 0.*alpha.*beta)./( 1-xK.*kappap-alpha.*delta_x ) );
index = delta_x==0 & residual_tes>=0 & residuals1<=0;
residuals1(index) = 0;

% (2) FOC of households on capital
residuals2 = muR + AL./p - rd - ( sigmaK + sigmap  ).^2.* yK - lambda.*( ( kappap - kappad ) ./ (  1-yK.*kappap-(1-yK).*kappad + kappa_fs )  ) ; 
residuals2( yK==0 ) = 0;

% (3) Equilibrium jump equations
w_post_jump = w .* (1-kappaw);
lambda_post_jump = lambda  + kappa_lambda_fun(lambda);
p_post_jump = model_sol.p_fitted( w_post_jump, lambda_post_jump );
residuals3 = p - p_post_jump - kappap;

% output
errors.residuals1=residuals1;
errors.residuals2=residuals2;
errors.residuals3=residuals3;


